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We briefly review a recently developed semiclassical theory 1 for quantum oscillations in 
the spatial (particle and kinetic energy) densities of finite fermion systems and present 
some examples of its results. We then discuss the inclusion of correlations (finite tem- 
peratures, pairing correlations) in the semiclassical theory. 

1. Introduction 

We have recently proposed 1 a semiclassical theory for quantum oscillations in the 
local particle densities and kinetic-energy densities of a system of N fermions in a 
local potential in D dimensions, described by the stationary Schrodinger equation 



The potential V(r) can be considered to represent the self-consistent mean field of 
an interacting system of fermions obtained in density functional theory (DFT). 2 
The single- particle wavefunctions n (r) then are the Kohn-Sham orbitals 3 and p(r) 
in (2) is the (ideally exact) ground-state particle density of the interacting system. 4 
Ordering the spectrum {E n } and choosing the energy scale such that < E\ < 
E% < ■ ■ ■ < E n < . . . , we fill the lowest levels up to the Fermi energy A and define 
the particle density by 



The factor 2 accounts for the spin degeneracy (the number N is assumed to be 
even). Further degeneracies, which may arise for systems in D > 1 dimensions, will 
not be spelled out bout included in the summations over n. For the kinetic-energy 
density, we consider two different definitions 
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E n <\ 





(3) 



E n <\ 



E n <\ 



which upon integration both yield the exact total kinetic energy. 
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The density of states g{E) of the system (1) is given by 

g(E)=Y^5{E-E n ), N = N(X)=2[ dEg(E). 

Jo 



(4) 



Separating its smooth and oscillatory parts by defining 

g(E):=g(E)+Sg(E), (5) 

the smooth part g(E) is given by the extended Thomas-Fermi (ETF) theory (see 
chapter 4.4.3 of Ref. 5 ), while the the oscillating part 8g(E) can be described, to 
leading order in h, by the semiclassical trace formula 6 ' 7 

1 



6g(E)~ ^ Apq{E) cos 



PO 



Spo(E) - — apo 



(6) 



The sum here is over all periodic orbits (POs) of the corresponding classical system 
described by the Hamilton function H(q, p) = p 2 /2m + V(q). Spo(E) is the action 
integral along the periodic orbit: 



Spo(E) = I p(£,q)-dq, 

JPO 



(7) 



with the classical momentum given by p(E, r) = (r / '|r|) ■J '2m[E — V^(r)]. For sys- 
tems in which all orbits are isolated in phase space, explicit expressions for the 
amplitudes Apo(E), which depend on the stabilities of the orbits, and for the 
Maslov indices crpo have been given by Gutzwiller. 6 For systems with continuous 
symmetries and for integrable systems, alternative expressions for the amplitudes 
and Maslov indices have been derived by many authors; they may be found in Ref. 5 
Separating smooth and oscillating terms of the spatial densities 

p(r) := p(r) + 6p(r) , r(r) := r (r) + <5r(r) , n (r) := ?i (r) + Sn (r) , (8) 

the smooth parts are given by the ETF theory. For their oscillating parts we have 
obtained 1 the following semiclassical expressions, valid again to leading order in h: 



5p(r) ~ ^ Ay(A,r) cos $ 7 (A, r) 



6t(t) 



6n(r) 



P 2 (A,r) 
2m 

P 2 (A,r) 
2m 



J^^l 7 (A, r) cos $ 7 (A,r) 



^^l 7 (A,r) Q 7 (A,r) cos $ 7 (A,r) 



(9) 
(10) 

(11) 



The sum here is over all closed orbits 7 starting and ending in the point r, and 

1 „ ,T n 7T 7T 



$ 7 (A,r) = -£ Y (A,r)--M 7 - i (D + l), 



(12) 



The action function S-y(X,r) = S-y(X,r,r' = r) is gained from the general open 
action integral for an orbit starting at r and ending at r' at fixed energy E = A 
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S 7 (A,r,r')=y p(A,q)-dq, (13) 

and /u 7 is the Morse index that counts the number of conjugate points along the 
orbit. 6 ' 7 For the functions _4 7 (A,r) and Q 7 (A, r) we refer to our articles. 1 ' 8 ' 9 The 
quantity A is the Fermi energy of the smooth (ETF) system, defined by 

X = X + 5X, TV = 2/ AEg{E). (14) 
Jo 

Since for POs the action integral Spo(A) is independent of r, they do not yield any 
oscillating phases in the above expressions; their contributions vary only smoothly 
with r through *4po(A,r) and Qpo(A, r). The leading contributions to the density 
oscillations come from the non-periodic orbits (NPOs). For one-dimensional sys- 
tems (D=l) it has, in fact, been shown 1 ' 8 that the contributions of the POs are 
completely absorbed by the smooth (TF) densities. In higher-dimensional systems, 
the POs must be included in (9) - (11) in connection with symmetry breaking at 
r = for spherical systems, and with bifurcations at finite distances |r| > in 
general, as demonstrated explicitly for the two-dimensional circular billiard. 9 



2. Selected results 

In this section we give some selected results of our semiclassical theory. We first 
present a very general result that may have interesting consequences for DFT. From 
(9), (10) one finds directly - without knowledge of the orbits 7 - the relation 

Sr(r) ~ [A - V(r)} 6p(r) , (15) 

which we call the (differential) local virial theorem (LVT) because it relates the 
potential and kinetic-energy densities locally at any given point r. The relation (15) 
was derived 10 for isotropic harmonic oscillators in arbitrary dimensions from their 
quantum- mechanical densities in the asymptotic limit TV — » 00. In our semiclassical 
theory it is obtained for arbitrary potentials. Since no assumption about the po- 
tential or the nature of the closed orbits 7 must be made to derive the LVT (15), 
it holds for arbitrary (integrable or non-integrable) systems in arbitrary dimensions 
with a local potential V(r), and hence also for interacting fermions in the mean- field 
approximation given by the DFT. We recall, however, that (15) is not expected to 
be valid close to the classical turning points where the semiclassical expressions (9) 
- (11) diverge and must be regularized by appropriate uniform approximations. 8,9 
A direct consequence of the LVT in (15) is the following relation: 

r(r) ~ r TF [p(r)]. (16) 

Hereby ttf [ptf (r)] = ttf(i") is the exact functional relation between the TF 
kinetic-energy and particle densities. Eq. (16) states that this TF functional (with- 
out gradient corrections!) holds approximately, for arbitrary local potentials V(r), 
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also between the exact quantum-mechanical densities r(r) and p(r) including their 
quantum oscillations. [It was shown in Ref. 8 to be exact up to first order in 6p(r).] 
In Fig. 1 wc test (16) explicitly for the coupled two-dimensional quartic oscillator 

V{x,y)= l -{x i + y i )-nx 2 y\ (17) 

whose classical dynamics is almost chaotic 11 ' 12 in the limits k = 1 and k — > — oo, 
but in practice also for k = 0.6 (see, e.g., Ref. 13 ). 
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Fig. 1. TF relation (16) for the potential (17) with re = 0.6 (units: h = m = 1) with N = 632. 
Cuts along the diagonal x = y. The solid line is the l.h.s., and the dashed line is the r.h.s. of (16). 



Wc find an excellent agreement over the whole region. That the TF kinetic-energy 
functional holds also for the oscillating exact densities to a surprising degree has 
been noted long ago, 14 but not understood until now. Similarly good numerical 
results are obtained also for the LVT (15), except very close to the classical turning 
points, for many systems 8 ' 9 with not too small particle numbers N. 




r/R 

Fig. 2. Particle density in the 2-dimensional circular billiard with various particle numbers N 
(units: h 2 /2m = R— 1). Dotted lines: quantum results, solid lines: semiclassical results. 9 
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Next we present some results for the particle densities. Figure 2 shows p(r) for 
four values of the number of TV particles bound in the two-dimensional circular 
billiard. The dotted line is the quantum result (2), and the solid line the converged 
semiclassical result (9), complemented by uniform approximations at the critical 
points as explained in detail in Ref. 9 Similar results are obtained also for the 
kinetic-energy densities, and for other types of potentials. 1,8 

It should be emphasized that, due to a factor 1/T 7 (A, r) in the semiclassical 
amplitudes Ay(A, r) in (9) - (11), the sums over the orbits 7 converge much faster 
than in the trace formula (6) for the level density. 

Note the Friedel oscillations in Fig. 2 near the surface (r = R), which are char- 
acteristic of a fermionic system near a steep boundary. In our semiclassical theory, 
the Friedel oscillations are caused by the shortest orbit with one reflection from the 
boundary (in Refs. 1,8 called the primitive "+" orbit). Its regularized contribution 
to the particle density of a spherical billiard in D dimensions is 8 

S P+ (r) = 2"I> + 1) J ^ , (18) 

where p^j is the TF density, v = D/2, z = 2 (R — r)p\/h, and p\ = (2mA) 1 / 2 is 
the smooth Fermi momentum. Integrating (18) over the whole space, we obtain 

5N + = J d D rSp + (r) = T^^" 1 ^ . ( 19 ) 

where So = 2n D / 2 R D ~ 1 /T(D/2) is the hypersurface of the Z?-dimensional sphere. 
It is interesting to note (see also Ref. 15 ) that (19) corresponds precisely to the 
surface term in the Weyl expansion 16 of the particle number N(X) which varies 
smoothly with the Fermi energy (the volume term being given by the TF theory). 

3. Inclusion of finite temperatures in the semiclassical theory 

In the following we outline how to include finite temperatures in the semiclassical 
formalism. Extensions of semiclassical trace formulae to finite temperatures have 
been used long ago in the context of nuclear physics 17 and more recently in meso- 
scopic physics. 18 We shall present here a derivation by means of a suitable folding 
function, which has proved useful also in the corresponding microscopic theories 19 
and allows for a straightforward generalization to include other types of correlations. 

For a grand-canonical ensemble of fcrmions embedded in a heat bath with fixed 
temperature, the variational energy is the so-called grand potential Q defined by 

n = (H)-TS- X(N) , (20) 

where H and N are the Hamilton and particle number operators, respectively, T is 
the temperature in energy units (i.e., we put the Boltzmann constant fce equal to 
unity), S is the entropy, and A the chemical potentials Note that both energy and 

a The quantities S and T without subscripts should not be confused with the actions Sj and 
periods T 7 of the classical orbits. 
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particle number are conserved only on the average. For non-interacting particles, 
we can write the Hclmholtz free energy F as 

F=(H)-TS = 2Y J E n v n -TS. (21) 

n 

Here v n are the Fermi occupation numbers 

^-XY . ( 22 ) 



1 + exp 

and the entropy S is given by 

S = -2 K log v n + {l- v n ) log(l -v n )]. (23) 

n 

The chemical potential A is determined by fixing the average particle number 

N = (N) = 2j2"n- (24) 

n 

Note that all sums in (21) - (24) and below run over the complete (infinite) spectrum 
of the Hamiltonian H . 

It has been shown 19 that the above quantities F, N and S can be expressed in 
terms of a convoluted finite-temperature level density gr{E) defined by a convolution 
of the "cold" (T = 0) density of states (4) 

/OO 
g(E') f T (E - E') dE' = ME - E n ) , (25) 
-oo „ 



whereby the folding function Jt{E) is given as 

1 

4T cosl?(£;/2T) 

The free energy then is given by 



ME) = T^ZTT^T^ ■ (26) 



F = 2 f Eg T (E)dE, (27) 

J — oo 

and the average particle number by 

N = 2 I 9t{E) dE . (28) 

J — oo 

To show that the integral (27) gives the correct free energy (21), including the "heat 
energy" —TS, requires some algebraic manipulations. From F, the entropy S can 
always be gained by the canonical relation 

dF 

S— w . (29) 

The same convolution can now be applied also to the semiclassical trace formula 
(6) for the oscillating part of the density of states which we re-write as 

Sg(E)~RcJ2A P o(E)e 1 *™^ (30) 

PO 
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with the phase 

$po(£) = ^Spo(£)-|<7p . (31) 

The oscillating part SgxiE) of the finite-temperature level density is obtained by 
the convolution of (30) with the function Jt{E) as in (25). In the spirit of the 
stationary-phase approximation, we take the slowly varying amplitude ^4po(£0 
outside of the integration and approximate the action in the phase by 

Spo(E') ~ Spo(E) + (E' - E) T PO {E) , (32) 

so that the result becomes a modified trace formula 

Sg T (E) ^RcJ2A P o(E)f T [ / rpo(E)}e^^ E \ (33) 

PO 

where 

T F o(E)=T PO (E)/h (34) 

and the temperature modulation factor fx is given by the Fourier transform of the 
convolution function fx'- 

/oc 
/r(w)e iTw dw. (35) 
-00 

The Fourier transform of the function (26) is known 20 and yields 

^ = ^ry (36 > 

The "hot" trace formulae (33) with the modulation factor (36) has previously been 
obtained in Refs. 17,18 The trace formula for the oscillating part of the free energy 
then becomes 5,17 to leading order in h 

6F~Re^A PO (X) (—^-) /t[T po (A)] (37) 
po V Tp °( A )/ 
For the spatial densities we can proceed exactly in the same way. For the particle 
density, e.g., the microscopic expression (2) is replaced by 

p T (r) = 2^|0„(r)| 2 ^, (38) 

n 

where the sum again runs over the complete spectrum. Starting from the semiclas- 
sical expression (9) for Sp(r) at T = 0, we rewrite it as 

5 Po (\,r) ~Re^^ 7 (A,r)e i *-'^' r) , (39) 

7 

where i' 7 (A, r) is the phase (12). The finite-T expression is given by the convolution 
integral 

Sp T (\,r)~ f S Po (X-E,r)f T {E)dE. (40) 

J —OO 
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Expanding the phase under the integral as above, we arrive at 

5p T (X,r) - Rc^^ 7 (A,r)/ T [T 7 (A,r)]e^^- r ) , (41) 

7 

where Xy(A,r) = T 7 (A, v)/h. The corresponding expressions for the temperature- 
dependent kinetic-energy densities are obvious. 

For the smooth parts of the densities, we recall that the ETF theory at T > 
is well known (see, e.g. Ref. 21 , where expressions up to 4-th order in h are given, 
and the literature quoted therein). 

Other types of correlations can be included in the semiclassical theory in the 
same way, as soon as a suitable folding function f CO rr{E) - corresponding to /t(E) 
in (26) - and its Fourier transform are known. One example is given by the pairing 
correlations discussed in the following section. 



4. Inclusion of pairing correlations in the BCS approximation 

A self-consistent microscopic approach to include pairing correlations is given by 
the Hartree-Fock-Bogolyubov (HFB) approach; we refer to an extended article 19 
for a recapitulation of this theory and the relevant literature. In the simplified BCS 
approach with constant paring gap A, the total energy of a system is written as 

-Bbcs = J2 E " v n - A J2 UnVn ' ( 42 ) 

n n 

where the sum goes over the complete spectrum (including all degeneracies) and 
the occupation numbers u n and v n are given by 

nV2 



Mi 



H (A-£„) 



VT^l. (43) 
Hereby £ n is the so-called quasiparticle energy 



£ n (A) = +V /(A - E n y + A 2 . (44) 

It was shown 19 that the BCS energy (42) is correctly given, including the pair 
condensation energy 

E p = - A ^ u nVn , (45) 
n 

by the convolution integral 

Ebcs = f Ef A (E)dE, (46) 

J — oo 

where the folding function /a(E) is defined as 

Ja(E) := A " -■„ ■ (47) 

2\E 2 + A 2 ] 3/2 
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The Fermi energy A in all above expressions is fixed by the average particle number: 

N= / g A (E)dE = J2<- (48) 

J-OD n 

The "paired" level density g A (E) is given by 

The Fourier transform of f A (E) is found 20 to be 

f A (T) = ATK^AT) , (50) 

where Ki(z) is a modified Bcsscl function. 22 Hence, replacing /t in (33) by /a, 
the semiclassical trace formula for the oscillating part of the paired level density 
becomes 

5g A (E) ~ReJ2Apo(E) fAlTpo(E)]e^ E \ (51) 

PO 

The trace formula for the oscillating part of the total BCS energy becomes, analo- 
gously to (37), 

8E BCS ~ Rc^po(A) ( — ^J] /a[T po (A)] e^^. (52) 
po \7po(A)/ 

That for the pair condensation energy, using E p = A dE^cs /dA and exploiting a 
recurrence relation for the Bessel functions, 22 becomes 

SE P ~A 2 ReJ2 -4po (A) K [T PO (A)] e 4$po & . (53) 

PO 

A similar result has recently been obtained in Ref. 23 

For the spatial densities we can, in principle, proceed as above. The pair- 
correlated particle density is quantum-mechanically given by 19 

PA (r)=J2\Ur)\ 2 v n . (54) 

n 

The semiclassical expression of its oscillating part becomes, similarly as above, 

6 PA (A, r) ~ Re £ A 1 (A, r) f A [T 7 (A)] . (55) 

7 

Corresponding results hold for the pair-correlated kinetic-energy densities. 

This is, however, not the end of the story. If one wants to express the pair- 
condensation energy (45) as a space integral, one requires an anomalous density 
matrix «(r,r'), defined by 19 

K (r,r') = ^0„(r)^(r'), (56) 

n 

where n refers to the time-reversed state of n. The semiclassical evaluation of this 
anomalous density matrix is the object of our ongoing research. 
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